# Check requisite packages are installed.
packages <- c(
"plotly",
"dplyr"
)
for (pkg in packages) {
library(pkg, character.only = TRUE)
}
What do the results look like?
dirViking <- file.path(
getwd(), "LCAB_LawMorton1996-NumericalPoolCommunityScaling"
)
dirVikingResults <- file.path(
dirViking, "results-2021-04"
)
resultFormat <- paste0(
"run-",
"%d", # Combination Number, or CombnNum.
"-",
"%s", # Run Seed.
".RDS"
)
# Copied from LawMorton1996-NumericalPoolCommunityScaling-Calculation.R
# TODO: In the future, make this a separate file that everyone can call...
set.seed(38427042)
basal <- c(3, 10, 30, 100, 300, 1000)
consumer <- c(3, 10, 30, 100, 300, 1000) * 2
events <- (max(basal) + max(consumer)) * 2
runs <- 100
logBodySize <- c(-2, -1, -1, 1) # Morton and Law 1997 version.
parameters <- c(0.01, 10, 0.5, 0.2, 100, 0.1)
# Need to rerun seedsPrep to get the random number generation right for seedsRun
seedsPrep <- runif(2 * length(basal) * length(consumer)) * 1E8
seedsRun <- runif(runs * length(basal) * length(consumer)) * 1E8
paramFrame <- with(list(
b = rep(basal, times = length(consumer)),
c = rep(consumer, each = length(basal)),
s1 = seedsPrep[1:(length(basal) * length(consumer))],
s2 = seedsPrep[
(length(basal) * length(consumer) + 1):(
2 * length(basal) * length(consumer))
],
sR = seedsRun
), {
temp <- data.frame(
CombnNum = 0,
Basals = b,
Consumers = c,
SeedPool = s1,
SeedMat = s2,
SeedRuns = "",
SeedRunsNum = 0,
EndStates = I(rep(list(""), length(b))),
EndStatesNum = 0,
EndStateSizes = I(rep(list(""), length(b))),
EndStateAssembly = I(rep(list(""), length(b)))
)
for (i in 1:nrow(temp)) {
seeds <- sR[((i - 1) * runs + 1) : (i * runs)]
temp$SeedRuns[i] <- toString(seeds) # CSV
temp$SeedRunsNum[i] <- length(seeds)
}
temp$CombnNum <- 1:nrow(temp)
temp
})
# Note: n + 2 end states. Failure to finish, failure to obtain state, and state.
for (i in 1:nrow(paramFrame)) {
resultsList <- list(
"No Run" = 0,
"No State" = 0
)
resultsSize <- list(
"0" = 0
)
resultsAssembly <- list(
"No Run" = data.frame(),
"No State" = data.frame()
)
seeds <- unlist(strsplit(paramFrame$SeedRuns[i], ', '))
for (seed in seeds) {
fileName <- file.path(
dirVikingResults,
sprintf(resultFormat, paramFrame$CombnNum[i], seed)
)
if (file.exists(fileName)) {
temp <- load(fileName)
temp <- eval(parse(text = temp)) # Get objects.
if (is.data.frame(temp)) {
community <- toString(
temp[[ncol(temp)]][[nrow(temp)]]
)
size <- toString(length(temp[[ncol(temp)]][[nrow(temp)]]))
if (community == "") {
resultsList$`No State` <- resultsList$`No State` + 1
resultsSize$`0` <- resultsSize$`0` + 1
} else if (community %in% names(resultsList)) {
resultsList[[community]] <- resultsList[[community]] + 1
resultsSize[[size]] <- resultsSize[[size]] + 1
} else {
resultsList[[community]] <- 1
resultsAssembly[[community]] <- temp
if (size %in% resultsSize) {
resultsSize[[size]] <- resultsSize[[size]] + 1
} else {
resultsSize[[size]] <- 1
}
}
} else {
resultsList$`No State` <- resultsList$`No State` + 1
resultsSize$`0` <- resultsSize$`0` + 1
}
} else {
resultsList$`No Run` <- resultsList$`No Run` + 1
resultsSize$`0` <- resultsSize$`0` + 1
}
}
paramFrame$EndStates[[i]] <- resultsList
paramFrame$EndStatesNum[i] <- length(resultsList) - 2 # ! No State, No Run
paramFrame$EndStateSizes[[i]] <- resultsSize
paramFrame$EndStateSizesNum[i] <- length(resultsSize) - 1 # ! 0
paramFrame$EndStateAssembly[[i]] <- resultsAssembly
}
# X, Y, Basal and Consumer.
# Z = Sizes of the Endstates.
plotScalingData <- data.frame(
CombnNum = rep(paramFrame$CombnNum, paramFrame$EndStatesNum),
Basals = rep(paramFrame$Basals, paramFrame$EndStatesNum),
Consumers = rep(paramFrame$Consumers, paramFrame$EndStatesNum)
)
# Communities
comms <- unlist(lapply(paramFrame$EndStates, names))
freqs <- unlist(paramFrame$EndStates)
asmbl <- unlist(paramFrame$EndStateAssembly, recursive = FALSE)
asmbl <- asmbl[comms != "No Run" & comms != "No State"]
freqs <- freqs[comms != "No Run" & comms != "No State"]
comms <- comms[comms != "No Run" & comms != "No State"]
asmbl <- lapply(asmbl, function(d) {
d %>% dplyr::filter(Result.Outcome != "Type 1 (Failure)" &
Result.Outcome != "Present")
})
plotScalingData$Communities <- comms
plotScalingData$CommunityFreq <- freqs
plotScalingData$CommunitySeq <- asmbl
# Community Size
temp <- unlist(lapply(strsplit(plotScalingData$Communities, ','), length))
plotScalingData$CommunitySize <- temp
# For usage by the reader.
plotScaling <- plotly::plot_ly(
plotScalingData,
x = ~Basals,
y = ~Consumers,
z = ~CommunitySize
)
plotScaling <- plotly::add_markers(plotScaling)
plotScaling <- plotly::layout(
plotScaling,
scene = list(
xaxis = list(type = "log"),
yaxis = list(type = "log")
)
)
plotScaling
plotScalingData %>% dplyr::select(-CommunitySeq)
How do they compare to each other?
# > runif(1) * 1E8
# [1] 82598679
set.seed(82598679)
temp <- load(file.path(
dirViking,
"LawMorton1996-NumericalPoolCommunityScaling-PoolMats.RDS"
))
mats <- eval(parse(text = temp[1]))
pools <- eval(parse(text = temp[2]))
candidateData <- plotScalingData %>% dplyr::group_by(
CombnNum
) %>% dplyr::mutate(
OtherSteadyStates = dplyr::n() - 1
) %>% dplyr::filter(
OtherSteadyStates > 0
)
candidateData %>% dplyr::select(-CommunitySeq)
candidateData$CommunityAbund <- ""
Warning messages:
1: Unknown or uninitialised column: `CommunityProd`.
2: Unknown or uninitialised column: `CommunityProd`.
3: Unknown or uninitialised column: `CommunityProd`.
for (r in 1:nrow(candidateData)) {
candidateData$CommunityAbund[r] <- toString(with(
candidateData[r, ], {
invasions <- CommunitySeq[[1]]$Result.Addition
abundance <- rep(0, Basals + Consumers)
for (i in invasions) {
abundance[i] <- abundance[i] + 1
abundance <- RMTRCode2::quiet(
deSolve::ode(
#rootSolve::steady(
y = abundance,
times = c(0:10000),
func = RMTRCode2::GeneralisedLotkaVolterra,
parms = list(a = mats[[CombnNum]],
r = pools[[CombnNum]]$ReproductionRate),
#method = "runsteady"
)[10001, -1])
}
#print(paste("Expected", Communities))
#print(paste("Calculated", toString(which(abundance > 0))))
#print(paste("> Epsilon:", toString(which(abundance > 1E-6))))
if (Communities == toString(which(abundance > 1E-6))) {
abundance[abundance > 1E-6]
} else {
"Failure"
}
}
))
}
candidateData$CommunityProd <- NA
for (r in 1:nrow(candidateData)) {
candidateData$CommunityProd[r] <- with(candidateData[r, ],
RMTRCode2::Productivity(
Pool = pools[[CombnNum]],
InteractionMatrix = mats[[CombnNum]],
Community = Communities,
Populations = CommunityAbund
)
)
}
islandFUN <- function(i, dat, pool, mat, dmat) {
temp <- dat[i, ]
RMTRCode2::IslandDynamics(
Pool = pool,
InteractionMatrix = mat,
Communities = c(
list(temp$Communities[1]),
rep("", nrow(dmat) - 2),
temp$Communities[2]
),
Populations = c(
list(temp$CommunityAbund[1]),
rep("", nrow(dmat) - 2),
list(temp$CommunityAbund[2])
),
DispersalPool = 0.0001,
DispersalIsland = dmat,
)
}
# For each group,
# For each pair,
# Run Island Dynamics,
# Save the result with its pairing
for (grp in unique(candidateData$CombnNum)) {
candidateDataSubset <- candidateData %>% dplyr::filter(CombnNum == grp)
pairingResults<- combn(
nrow(candidateDataSubset), 2,
islandFUN,
dat = candidateDataSubset,
pool = pools[[grp]],
mat = mats[[grp]],
dmat = matrix(c(0, 1, 1, 0), nrow = 2, ncol = 2),
simplify = FALSE
)
}
for (grp in unique(candidateData$CombnNum)) {
candidateDataSubset <- candidateData %>% dplyr::filter(CombnNum == grp)
pairingResults<- combn(
nrow(candidateDataSubset), 2,
islandFUN,
dat = candidateDataSubset,
pool = pools[[grp]],
mat = mats[[grp]],
dmat = matrix(c(
0, 1, 0, # Island 2 -> 1
1, 0, 1, # Island 1 -> 2, Island 3 -> 2
0, 1, 0 # Island 2 -> 3
), nrow = 3, ncol = 3, byrow = TRUE),
simplify = FALSE
)
}
LS0tDQp0aXRsZTogIlZpa2luZyBSZXN1bHRzLCAyMDIxLTA0Ig0Kb3V0cHV0Og0KICBodG1sX25vdGVib29rOg0KICAgIGNvZGVfZm9sZGluZzogaGlkZQ0KLS0tDQoNCmBgYHtyIGxpYnN9DQojIENoZWNrIHJlcXVpc2l0ZSBwYWNrYWdlcyBhcmUgaW5zdGFsbGVkLg0KcGFja2FnZXMgPC0gYygNCiAgInBsb3RseSIsIA0KICAiZHBseXIiDQopDQpmb3IgKHBrZyBpbiBwYWNrYWdlcykgew0KICBsaWJyYXJ5KHBrZywgY2hhcmFjdGVyLm9ubHkgPSBUUlVFKQ0KfQ0KYGBgDQoNCiMgV2hhdCBkbyB0aGUgcmVzdWx0cyBsb29rIGxpa2U/DQpgYGB7ciBkaXJzfQ0KZGlyVmlraW5nIDwtIGZpbGUucGF0aCgNCiAgZ2V0d2QoKSwgIkxDQUJfTGF3TW9ydG9uMTk5Ni1OdW1lcmljYWxQb29sQ29tbXVuaXR5U2NhbGluZyINCikNCmRpclZpa2luZ1Jlc3VsdHMgPC0gZmlsZS5wYXRoKA0KICBkaXJWaWtpbmcsICJyZXN1bHRzLTIwMjEtMDQiDQopDQpyZXN1bHRGb3JtYXQgPC0gcGFzdGUwKA0KICAicnVuLSIsIA0KICAiJWQiLCAjIENvbWJpbmF0aW9uIE51bWJlciwgb3IgQ29tYm5OdW0uDQogICItIiwgDQogICIlcyIsICMgUnVuIFNlZWQuDQogICIuUkRTIg0KKQ0KYGBgDQoNCmBgYHtyIHBhcmFtc30NCiMgQ29waWVkIGZyb20gTGF3TW9ydG9uMTk5Ni1OdW1lcmljYWxQb29sQ29tbXVuaXR5U2NhbGluZy1DYWxjdWxhdGlvbi5SDQojIFRPRE86IEluIHRoZSBmdXR1cmUsIG1ha2UgdGhpcyBhIHNlcGFyYXRlIGZpbGUgdGhhdCBldmVyeW9uZSBjYW4gY2FsbC4uLg0Kc2V0LnNlZWQoMzg0MjcwNDIpDQoNCmJhc2FsIDwtIGMoMywgMTAsIDMwLCAxMDAsIDMwMCwgMTAwMCkNCmNvbnN1bWVyIDwtIGMoMywgMTAsIDMwLCAxMDAsIDMwMCwgMTAwMCkgKiAyDQpldmVudHMgPC0gKG1heChiYXNhbCkgKyBtYXgoY29uc3VtZXIpKSAqIDINCnJ1bnMgPC0gMTAwDQoNCmxvZ0JvZHlTaXplIDwtIGMoLTIsIC0xLCAtMSwgMSkgIyBNb3J0b24gYW5kIExhdyAxOTk3IHZlcnNpb24uDQpwYXJhbWV0ZXJzIDwtIGMoMC4wMSwgMTAsIDAuNSwgMC4yLCAxMDAsIDAuMSkNCg0KIyBOZWVkIHRvIHJlcnVuIHNlZWRzUHJlcCB0byBnZXQgdGhlIHJhbmRvbSBudW1iZXIgZ2VuZXJhdGlvbiByaWdodCBmb3Igc2VlZHNSdW4NCnNlZWRzUHJlcCA8LSBydW5pZigyICogbGVuZ3RoKGJhc2FsKSAqIGxlbmd0aChjb25zdW1lcikpICogMUU4DQpzZWVkc1J1biA8LSBydW5pZihydW5zICogbGVuZ3RoKGJhc2FsKSAqIGxlbmd0aChjb25zdW1lcikpICogMUU4DQpgYGANCg0KYGBge3Igb3JnYW5pc2VQYXJhbXN9DQpwYXJhbUZyYW1lIDwtIHdpdGgobGlzdCgNCiAgYiA9IHJlcChiYXNhbCwgdGltZXMgPSBsZW5ndGgoY29uc3VtZXIpKSwNCiAgYyA9IHJlcChjb25zdW1lciwgZWFjaCA9IGxlbmd0aChiYXNhbCkpLA0KICBzMSA9IHNlZWRzUHJlcFsxOihsZW5ndGgoYmFzYWwpICogbGVuZ3RoKGNvbnN1bWVyKSldLA0KICBzMiA9IHNlZWRzUHJlcFsNCiAgICAobGVuZ3RoKGJhc2FsKSAqIGxlbmd0aChjb25zdW1lcikgKyAxKTooDQogICAgICAyICogbGVuZ3RoKGJhc2FsKSAqIGxlbmd0aChjb25zdW1lcikpDQogIF0sDQogIHNSID0gc2VlZHNSdW4NCiksIHsNCiAgdGVtcCA8LSBkYXRhLmZyYW1lKA0KICAgIENvbWJuTnVtID0gMCwNCiAgICBCYXNhbHMgPSBiLA0KICAgIENvbnN1bWVycyA9IGMsDQogICAgU2VlZFBvb2wgPSBzMSwNCiAgICBTZWVkTWF0ID0gczIsDQogICAgU2VlZFJ1bnMgPSAiIiwNCiAgICBTZWVkUnVuc051bSA9IDAsDQogICAgRW5kU3RhdGVzID0gSShyZXAobGlzdCgiIiksIGxlbmd0aChiKSkpLA0KICAgIEVuZFN0YXRlc051bSA9IDAsDQogICAgRW5kU3RhdGVTaXplcyA9IEkocmVwKGxpc3QoIiIpLCBsZW5ndGgoYikpKSwNCiAgICBFbmRTdGF0ZUFzc2VtYmx5ID0gSShyZXAobGlzdCgiIiksIGxlbmd0aChiKSkpDQogICkNCiAgZm9yIChpIGluIDE6bnJvdyh0ZW1wKSkgew0KICAgIHNlZWRzIDwtIHNSWygoaSAtIDEpICogcnVucyArIDEpIDogKGkgKiBydW5zKV0NCiAgICB0ZW1wJFNlZWRSdW5zW2ldIDwtIHRvU3RyaW5nKHNlZWRzKSAjIENTVg0KICAgIHRlbXAkU2VlZFJ1bnNOdW1baV0gPC0gbGVuZ3RoKHNlZWRzKQ0KICB9DQogIHRlbXAkQ29tYm5OdW0gPC0gMTpucm93KHRlbXApDQogIHRlbXANCn0pDQpgYGANCg0KYGBge3IgbG9hZFJlc3VsdHN9DQojIE5vdGU6IG4gKyAyIGVuZCBzdGF0ZXMuIEZhaWx1cmUgdG8gZmluaXNoLCBmYWlsdXJlIHRvIG9idGFpbiBzdGF0ZSwgYW5kIHN0YXRlLg0KZm9yIChpIGluIDE6bnJvdyhwYXJhbUZyYW1lKSkgew0KICByZXN1bHRzTGlzdCA8LSBsaXN0KA0KICAgICJObyBSdW4iID0gMCwNCiAgICAiTm8gU3RhdGUiID0gMA0KICApDQogIHJlc3VsdHNTaXplIDwtIGxpc3QoDQogICAgIjAiID0gMA0KICApDQogIHJlc3VsdHNBc3NlbWJseSA8LSBsaXN0KA0KICAgICJObyBSdW4iID0gZGF0YS5mcmFtZSgpLA0KICAgICJObyBTdGF0ZSIgPSBkYXRhLmZyYW1lKCkNCiAgKQ0KICBzZWVkcyA8LSB1bmxpc3Qoc3Ryc3BsaXQocGFyYW1GcmFtZSRTZWVkUnVuc1tpXSwgJywgJykpDQogIGZvciAoc2VlZCBpbiBzZWVkcykgew0KICAgIGZpbGVOYW1lIDwtIGZpbGUucGF0aCgNCiAgICAgIGRpclZpa2luZ1Jlc3VsdHMsDQogICAgICBzcHJpbnRmKHJlc3VsdEZvcm1hdCwgcGFyYW1GcmFtZSRDb21ibk51bVtpXSwgc2VlZCkNCiAgICApDQogICAgDQogICAgaWYgKGZpbGUuZXhpc3RzKGZpbGVOYW1lKSkgew0KICAgICAgdGVtcCA8LSBsb2FkKGZpbGVOYW1lKQ0KICAgICAgdGVtcCA8LSBldmFsKHBhcnNlKHRleHQgPSB0ZW1wKSkgIyBHZXQgb2JqZWN0cy4NCiAgICAgIA0KICAgICAgaWYgKGlzLmRhdGEuZnJhbWUodGVtcCkpIHsNCiAgICAgICAgY29tbXVuaXR5IDwtIHRvU3RyaW5nKA0KICAgICAgICAgIHRlbXBbW25jb2wodGVtcCldXVtbbnJvdyh0ZW1wKV1dDQogICAgICAgICkNCiAgICAgICAgc2l6ZSA8LSB0b1N0cmluZyhsZW5ndGgodGVtcFtbbmNvbCh0ZW1wKV1dW1tucm93KHRlbXApXV0pKQ0KICAgICAgICANCiAgICAgICAgaWYgKGNvbW11bml0eSA9PSAiIikgew0KICAgICAgICAgIHJlc3VsdHNMaXN0JGBObyBTdGF0ZWAgPC0gcmVzdWx0c0xpc3QkYE5vIFN0YXRlYCArIDENCiAgICAgICAgICByZXN1bHRzU2l6ZSRgMGAgPC0gcmVzdWx0c1NpemUkYDBgICsgMQ0KICAgICAgICAgIA0KICAgICAgICB9IGVsc2UgaWYgKGNvbW11bml0eSAlaW4lIG5hbWVzKHJlc3VsdHNMaXN0KSkgew0KICAgICAgICAgIHJlc3VsdHNMaXN0W1tjb21tdW5pdHldXSA8LSByZXN1bHRzTGlzdFtbY29tbXVuaXR5XV0gKyAxDQogICAgICAgICAgcmVzdWx0c1NpemVbW3NpemVdXSA8LSByZXN1bHRzU2l6ZVtbc2l6ZV1dICsgMQ0KICAgICAgICAgIA0KICAgICAgICB9IGVsc2Ugew0KICAgICAgICAgIHJlc3VsdHNMaXN0W1tjb21tdW5pdHldXSA8LSAxDQogICAgICAgICAgcmVzdWx0c0Fzc2VtYmx5W1tjb21tdW5pdHldXSA8LSB0ZW1wDQogICAgICAgICAgDQogICAgICAgICAgaWYgKHNpemUgJWluJSByZXN1bHRzU2l6ZSkgew0KICAgICAgICAgICAgcmVzdWx0c1NpemVbW3NpemVdXSA8LSByZXN1bHRzU2l6ZVtbc2l6ZV1dICsgMQ0KICAgICAgICAgIH0gZWxzZSB7DQogICAgICAgICAgICByZXN1bHRzU2l6ZVtbc2l6ZV1dIDwtIDENCiAgICAgICAgICB9DQogICAgICAgIH0NCiAgICAgIH0gZWxzZSB7DQogICAgICAgIHJlc3VsdHNMaXN0JGBObyBTdGF0ZWAgPC0gcmVzdWx0c0xpc3QkYE5vIFN0YXRlYCArIDENCiAgICAgICAgcmVzdWx0c1NpemUkYDBgIDwtIHJlc3VsdHNTaXplJGAwYCArIDENCiAgICAgIH0NCiAgICB9IGVsc2Ugew0KICAgICAgcmVzdWx0c0xpc3QkYE5vIFJ1bmAgPC0gcmVzdWx0c0xpc3QkYE5vIFJ1bmAgKyAxDQogICAgICByZXN1bHRzU2l6ZSRgMGAgPC0gcmVzdWx0c1NpemUkYDBgICsgMQ0KICAgIH0NCiAgfQ0KICANCiAgcGFyYW1GcmFtZSRFbmRTdGF0ZXNbW2ldXSA8LSByZXN1bHRzTGlzdA0KICBwYXJhbUZyYW1lJEVuZFN0YXRlc051bVtpXSA8LSBsZW5ndGgocmVzdWx0c0xpc3QpIC0gMiAjICEgTm8gU3RhdGUsIE5vIFJ1bg0KICBwYXJhbUZyYW1lJEVuZFN0YXRlU2l6ZXNbW2ldXSA8LSByZXN1bHRzU2l6ZQ0KICBwYXJhbUZyYW1lJEVuZFN0YXRlU2l6ZXNOdW1baV0gPC0gbGVuZ3RoKHJlc3VsdHNTaXplKSAtIDEgIyAhIDANCiAgcGFyYW1GcmFtZSRFbmRTdGF0ZUFzc2VtYmx5W1tpXV0gPC0gcmVzdWx0c0Fzc2VtYmx5DQp9DQpgYGANCg0KPCEtLWBgYHtyIHNob3dSZXN1bHRzfQ0KcGFyYW1GcmFtZVssIGMoMTozLCA4OjEwKV0NCmBgYC0tPg0KDQpgYGB7ciBwbG90M0R9DQojIFgsIFksIEJhc2FsIGFuZCBDb25zdW1lci4NCiMgWiA9IFNpemVzIG9mIHRoZSBFbmRzdGF0ZXMuDQoNCnBsb3RTY2FsaW5nRGF0YSA8LSBkYXRhLmZyYW1lKA0KICBDb21ibk51bSA9IHJlcChwYXJhbUZyYW1lJENvbWJuTnVtLCBwYXJhbUZyYW1lJEVuZFN0YXRlc051bSksDQogIEJhc2FscyA9IHJlcChwYXJhbUZyYW1lJEJhc2FscywgcGFyYW1GcmFtZSRFbmRTdGF0ZXNOdW0pLA0KICBDb25zdW1lcnMgPSByZXAocGFyYW1GcmFtZSRDb25zdW1lcnMsIHBhcmFtRnJhbWUkRW5kU3RhdGVzTnVtKQ0KKQ0KDQojIENvbW11bml0aWVzDQpjb21tcyA8LSB1bmxpc3QobGFwcGx5KHBhcmFtRnJhbWUkRW5kU3RhdGVzLCBuYW1lcykpDQpmcmVxcyA8LSB1bmxpc3QocGFyYW1GcmFtZSRFbmRTdGF0ZXMpDQphc21ibCA8LSB1bmxpc3QocGFyYW1GcmFtZSRFbmRTdGF0ZUFzc2VtYmx5LCByZWN1cnNpdmUgPSBGQUxTRSkNCmFzbWJsIDwtIGFzbWJsW2NvbW1zICE9ICJObyBSdW4iICYgY29tbXMgIT0gIk5vIFN0YXRlIl0NCmZyZXFzIDwtIGZyZXFzW2NvbW1zICE9ICJObyBSdW4iICYgY29tbXMgIT0gIk5vIFN0YXRlIl0NCmNvbW1zIDwtIGNvbW1zW2NvbW1zICE9ICJObyBSdW4iICYgY29tbXMgIT0gIk5vIFN0YXRlIl0NCg0KYXNtYmwgPC0gbGFwcGx5KGFzbWJsLCBmdW5jdGlvbihkKSB7DQogIGQgJT4lIGRwbHlyOjpmaWx0ZXIoUmVzdWx0Lk91dGNvbWUgIT0gIlR5cGUgMSAoRmFpbHVyZSkiICYgDQogICAgICAgICAgICAgICAgICAgICAgICBSZXN1bHQuT3V0Y29tZSAhPSAiUHJlc2VudCIpDQp9KQ0KDQpwbG90U2NhbGluZ0RhdGEkQ29tbXVuaXRpZXMgPC0gY29tbXMNCnBsb3RTY2FsaW5nRGF0YSRDb21tdW5pdHlGcmVxIDwtIGZyZXFzDQpwbG90U2NhbGluZ0RhdGEkQ29tbXVuaXR5U2VxIDwtIGFzbWJsDQoNCiMgQ29tbXVuaXR5IFNpemUNCnRlbXAgPC0gdW5saXN0KGxhcHBseShzdHJzcGxpdChwbG90U2NhbGluZ0RhdGEkQ29tbXVuaXRpZXMsICcsJyksIGxlbmd0aCkpDQpwbG90U2NhbGluZ0RhdGEkQ29tbXVuaXR5U2l6ZSA8LSB0ZW1wDQoNCiMgRm9yIHVzYWdlIGJ5IHRoZSByZWFkZXIuDQoNCnBsb3RTY2FsaW5nIDwtIHBsb3RseTo6cGxvdF9seSgNCiAgcGxvdFNjYWxpbmdEYXRhLA0KICB4ID0gfkJhc2FscywNCiAgeSA9IH5Db25zdW1lcnMsDQogIHogPSB+Q29tbXVuaXR5U2l6ZQ0KKQ0KDQpwbG90U2NhbGluZyA8LSBwbG90bHk6OmFkZF9tYXJrZXJzKHBsb3RTY2FsaW5nKQ0KDQpwbG90U2NhbGluZyA8LSBwbG90bHk6OmxheW91dCgNCiAgcGxvdFNjYWxpbmcsDQogIHNjZW5lID0gbGlzdCgNCiAgICB4YXhpcyA9IGxpc3QodHlwZSA9ICJsb2ciKSwNCiAgICB5YXhpcyA9IGxpc3QodHlwZSA9ICJsb2ciKQ0KICApDQopDQoNCnBsb3RTY2FsaW5nDQpgYGANCmBgYHtyIHBsb3QzZERhdGF9DQpwbG90U2NhbGluZ0RhdGEgJT4lIGRwbHlyOjpzZWxlY3QoLUNvbW11bml0eVNlcSkNCmBgYA0KIyBIb3cgZG8gdGhleSBjb21wYXJlIHRvIGVhY2ggb3RoZXI/DQoNCmBgYHtyIGxvYWRQb29sc01hdHN9DQojID4gcnVuaWYoMSkgKiAxRTgNCiMgWzFdIDgyNTk4Njc5DQpzZXQuc2VlZCg4MjU5ODY3OSkNCg0KdGVtcCA8LSBsb2FkKGZpbGUucGF0aCgNCiAgZGlyVmlraW5nLCANCiAgIkxhd01vcnRvbjE5OTYtTnVtZXJpY2FsUG9vbENvbW11bml0eVNjYWxpbmctUG9vbE1hdHMuUkRTIg0KKSkNCm1hdHMgIDwtIGV2YWwocGFyc2UodGV4dCA9IHRlbXBbMV0pKQ0KcG9vbHMgPC0gZXZhbChwYXJzZSh0ZXh0ID0gdGVtcFsyXSkpDQpgYGANCg0KYGBge3IgY29tcHV0ZUNhbmRpZGF0ZXN9DQpjYW5kaWRhdGVEYXRhIDwtIHBsb3RTY2FsaW5nRGF0YSAlPiUgZHBseXI6Omdyb3VwX2J5KA0KICBDb21ibk51bQ0KKSAlPiUgZHBseXI6Om11dGF0ZSgNCiAgT3RoZXJTdGVhZHlTdGF0ZXMgPSBkcGx5cjo6bigpIC0gMQ0KKSAlPiUgZHBseXI6OmZpbHRlcigNCiAgT3RoZXJTdGVhZHlTdGF0ZXMgPiAwDQopDQpjYW5kaWRhdGVEYXRhICU+JSBkcGx5cjo6c2VsZWN0KC1Db21tdW5pdHlTZXEpDQpgYGANCg0KPCEtLWBgYHtyIGNvbXB1dGVQb3B1bGF0aW9uc30NCnBsb3RTY2FsaW5nRGF0YSRDb21tdW5pdHlBYnVuZCA8LSAiIg0KZm9yIChyIGluIDE6bnJvdyhwbG90U2NhbGluZ0RhdGEpKSB7DQogIHByaW50KHIpDQogIHBsb3RTY2FsaW5nRGF0YSRDb21tdW5pdHlBYnVuZFtyXSA8LSB3aXRoKHBsb3RTY2FsaW5nRGF0YVtyLCBdLCB0b1N0cmluZygNCiAgICBSTVRSQ29kZTI6OkZpbmRTdGVhZHlTdGF0ZUZyb21TaXplKA0KICAgICAgUG9vbCA9IHBvb2xzW1tDb21ibk51bV1dLCANCiAgICAgIEludGVyYWN0aW9uTWF0cml4ID0gbWF0c1tbQ29tYm5OdW1dXSwgDQogICAgICBDb21tdW5pdHkgPSBDb21tdW5pdGllcywgDQogICAgICBQb3B1bGF0aW9ucyA9IHJlcCgxMDAsIENvbW11bml0eVNpemUpLCAjIE5vIGdvb2QgZ3Vlc3NlcyBoZXJlDQogICAgICBtYXhSYW5kVmFsID0gMkU0LCAjIDIgKiByb3VuZCBvZiB0aGUgbGFyZ2VzdCB2YWx1ZSB3ZSBoYXZlIHNlZW4gc28gZmFyLg0KICAgICAgTWF4QXR0ZW1wdHMgPSAxRTQsDQogICAgICBWZXJib3NlID0gRkFMU0UNCiAgICApDQogICkpDQp9DQpgYGAtLT4NCg0KYGBge3IgY29tcHV0ZVBvcHVsYXRpb25zMn0NCmNhbmRpZGF0ZURhdGEkQ29tbXVuaXR5QWJ1bmQgPC0gIiINCmZvciAociBpbiAxOm5yb3coY2FuZGlkYXRlRGF0YSkpIHsNCiAgY2FuZGlkYXRlRGF0YSRDb21tdW5pdHlBYnVuZFtyXSA8LSB0b1N0cmluZyh3aXRoKA0KICAgIGNhbmRpZGF0ZURhdGFbciwgXSwgew0KICAgICAgaW52YXNpb25zIDwtIENvbW11bml0eVNlcVtbMV1dJFJlc3VsdC5BZGRpdGlvbg0KICAgICAgYWJ1bmRhbmNlIDwtIHJlcCgwLCBCYXNhbHMgKyBDb25zdW1lcnMpDQogICAgICBmb3IgKGkgaW4gaW52YXNpb25zKSB7DQogICAgICAgIGFidW5kYW5jZVtpXSA8LSBhYnVuZGFuY2VbaV0gKyAxDQogICAgICAgIGFidW5kYW5jZSA8LSBSTVRSQ29kZTI6OnF1aWV0KA0KICAgICAgICAgIGRlU29sdmU6Om9kZSgNCiAgICAgICAgICAgICNyb290U29sdmU6OnN0ZWFkeSgNCiAgICAgICAgICAgIHkgPSBhYnVuZGFuY2UsDQogICAgICAgICAgICB0aW1lcyA9IGMoMDoxMDAwMCksDQogICAgICAgICAgICBmdW5jID0gUk1UUkNvZGUyOjpHZW5lcmFsaXNlZExvdGthVm9sdGVycmEsDQogICAgICAgICAgICBwYXJtcyA9IGxpc3QoYSA9IG1hdHNbW0NvbWJuTnVtXV0sDQogICAgICAgICAgICAgICAgICAgICAgICAgciA9IHBvb2xzW1tDb21ibk51bV1dJFJlcHJvZHVjdGlvblJhdGUpLA0KICAgICAgICAgICAgI21ldGhvZCA9ICJydW5zdGVhZHkiDQogICAgICAgICAgKVsxMDAwMSwgLTFdKQ0KICAgICAgfQ0KICAgICAgI3ByaW50KHBhc3RlKCJFeHBlY3RlZCIsIENvbW11bml0aWVzKSkNCiAgICAgICNwcmludChwYXN0ZSgiQ2FsY3VsYXRlZCIsIHRvU3RyaW5nKHdoaWNoKGFidW5kYW5jZSA+IDApKSkpDQogICAgICAjcHJpbnQocGFzdGUoIj4gRXBzaWxvbjoiLCB0b1N0cmluZyh3aGljaChhYnVuZGFuY2UgPiAxRS02KSkpKQ0KICAgICAgaWYgKENvbW11bml0aWVzID09IHRvU3RyaW5nKHdoaWNoKGFidW5kYW5jZSA+IDFFLTYpKSkgew0KICAgICAgICBhYnVuZGFuY2VbYWJ1bmRhbmNlID4gMUUtNl0NCiAgICAgIH0gZWxzZSB7DQogICAgICAgICJGYWlsdXJlIg0KICAgICAgfQ0KICAgIH0NCiAgKSkNCn0NCg0KYGBgDQoNCmBgYHtyIGNvbXB1dGVQcm9kdWN0aXZpdHl9DQpjYW5kaWRhdGVEYXRhJENvbW11bml0eVByb2QgPC0gTkENCmZvciAociBpbiAxOm5yb3coY2FuZGlkYXRlRGF0YSkpIHsNCiAgY2FuZGlkYXRlRGF0YSRDb21tdW5pdHlQcm9kW3JdIDwtIHdpdGgoY2FuZGlkYXRlRGF0YVtyLCBdLCANCiAgICBSTVRSQ29kZTI6OlByb2R1Y3Rpdml0eSgNCiAgICAgIFBvb2wgPSBwb29sc1tbQ29tYm5OdW1dXSwgDQogICAgICBJbnRlcmFjdGlvbk1hdHJpeCA9IG1hdHNbW0NvbWJuTnVtXV0sIA0KICAgICAgQ29tbXVuaXR5ID0gQ29tbXVuaXRpZXMsIA0KICAgICAgUG9wdWxhdGlvbnMgPSBDb21tdW5pdHlBYnVuZA0KICAgICkNCiAgKQ0KfQ0KYGBgDQoNCmBgYHtyIGlzbGFuZEZVTn0NCmlzbGFuZEZVTiA8LSBmdW5jdGlvbihpLCBkYXQsIHBvb2wsIG1hdCwgZG1hdCkgew0KICAgICAgdGVtcCA8LSBkYXRbaSwgXQ0KICAgICAgUk1UUkNvZGUyOjpJc2xhbmREeW5hbWljcygNCiAgICAgICAgUG9vbCA9IHBvb2wsDQogICAgICAgIEludGVyYWN0aW9uTWF0cml4ID0gbWF0LA0KICAgICAgICBDb21tdW5pdGllcyA9IGMoDQogICAgICAgICAgbGlzdCh0ZW1wJENvbW11bml0aWVzWzFdKSwNCiAgICAgICAgICByZXAoIiIsIG5yb3coZG1hdCkgLSAyKSwNCiAgICAgICAgICB0ZW1wJENvbW11bml0aWVzWzJdDQogICAgICAgICksDQogICAgICAgIFBvcHVsYXRpb25zID0gYygNCiAgICAgICAgICBsaXN0KHRlbXAkQ29tbXVuaXR5QWJ1bmRbMV0pLA0KICAgICAgICAgIHJlcCgiIiwgbnJvdyhkbWF0KSAtIDIpLA0KICAgICAgICAgIGxpc3QodGVtcCRDb21tdW5pdHlBYnVuZFsyXSkNCiAgICAgICAgKSwNCiAgICAgICAgRGlzcGVyc2FsUG9vbCA9IDAuMDAwMSwNCiAgICAgICAgRGlzcGVyc2FsSXNsYW5kID0gZG1hdCwNCiAgICAgICkNCiAgICB9DQpgYGANCg0KYGBge3IgaXNsYW5kT25lVHdvfQ0KIyBGb3IgZWFjaCBncm91cCwNCiMgRm9yIGVhY2ggcGFpciwNCiMgUnVuIElzbGFuZCBEeW5hbWljcywNCiMgU2F2ZSB0aGUgcmVzdWx0IHdpdGggaXRzIHBhaXJpbmcNCg0KZm9yIChncnAgaW4gdW5pcXVlKGNhbmRpZGF0ZURhdGEkQ29tYm5OdW0pKSB7DQogIGNhbmRpZGF0ZURhdGFTdWJzZXQgPC0gY2FuZGlkYXRlRGF0YSAlPiUgZHBseXI6OmZpbHRlcihDb21ibk51bSA9PSBncnApDQogIHBhaXJpbmdSZXN1bHRzPC0gY29tYm4oDQogICAgbnJvdyhjYW5kaWRhdGVEYXRhU3Vic2V0KSwgMiwgDQogICAgaXNsYW5kRlVOLA0KICAgIGRhdCA9IGNhbmRpZGF0ZURhdGFTdWJzZXQsIA0KICAgIHBvb2wgPSBwb29sc1tbZ3JwXV0sDQogICAgbWF0ID0gbWF0c1tbZ3JwXV0sDQogICAgZG1hdCA9IG1hdHJpeChjKDAsIDEsIDEsIDApLCBucm93ID0gMiwgbmNvbCA9IDIpLA0KICAgIHNpbXBsaWZ5ID0gRkFMU0UNCiAgKQ0KfQ0KYGBgDQoNCmBgYHtyIGlzbGFuZE9uZUVtcHR5VHdvfQ0KZm9yIChncnAgaW4gdW5pcXVlKGNhbmRpZGF0ZURhdGEkQ29tYm5OdW0pKSB7DQogIGNhbmRpZGF0ZURhdGFTdWJzZXQgPC0gY2FuZGlkYXRlRGF0YSAlPiUgZHBseXI6OmZpbHRlcihDb21ibk51bSA9PSBncnApDQogIHBhaXJpbmdSZXN1bHRzPC0gY29tYm4oDQogICAgbnJvdyhjYW5kaWRhdGVEYXRhU3Vic2V0KSwgMiwgDQogICAgaXNsYW5kRlVOLA0KICAgIGRhdCA9IGNhbmRpZGF0ZURhdGFTdWJzZXQsIA0KICAgIHBvb2wgPSBwb29sc1tbZ3JwXV0sDQogICAgbWF0ID0gbWF0c1tbZ3JwXV0sDQogICAgZG1hdCA9IG1hdHJpeChjKA0KICAgICAgMCwgMSwgMCwgIyBJc2xhbmQgMiAtPiAxDQogICAgICAxLCAwLCAxLCAjIElzbGFuZCAxIC0+IDIsIElzbGFuZCAzIC0+IDINCiAgICAgIDAsIDEsIDAgICMgSXNsYW5kIDIgLT4gMw0KICAgICksIG5yb3cgPSAzLCBuY29sID0gMywgYnlyb3cgPSBUUlVFKSwNCiAgICBzaW1wbGlmeSA9IEZBTFNFDQogICkNCn0NCmBgYA0K